%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% MTO Exercises
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear
clc

% Load parameters
load param

x_s = myx;

% Runt MTO counterfcatuals
data_set5 = counter_fact_mto(x_s, [], 0.05,0);   % 5% scale
data_set10 = counter_fact_mto(x_s, [], 0.10,0);  % 10% scale
data_set15 = counter_fact_mto(x_s, [], 0.15,0);  % 15% scale
data_set20 = counter_fact_mto(x_s, [], 0.20,0);  % 20% scale
data_set25 = counter_fact_mto(x_s, [], 0.25,0);  % 25% scale
data_set0 = counter_fact_mto(x_s, [], 0.0015,0); % MTO scale

data_set_base = counter_fact_mto(x_s, [], 0,0); % 0%

save data_mto data_set_base data_set5 data_set10 data_set15 data_set20 data_set25  data_set0

Ss_ss = cell2mat(data_set_base(:,3))';
save Ss_ss Ss_ss
rs_ss = cell2mat(data_set_base(:,5))';
save rs_ss rs_ss
tax_ss = [cell2mat(data_set0(:,11)) cell2mat(data_set5(:,11)) cell2mat(data_set10(:,11)) cell2mat(data_set15(:,11)) cell2mat(data_set20(:,11)) cell2mat(data_set25(:,11))]';
save tax_ss tax_ss

vouch_share = [cell2mat(data_set0(1:5,13)) cell2mat(data_set5(1:5,13)) cell2mat(data_set10(1:5,13)) cell2mat(data_set15(1:5,13)) cell2mat(data_set20(1:5,13)) cell2mat(data_set25(1:5,13))]
takeup      = [cell2mat(data_set0(:,14)) cell2mat(data_set5(:,14)) cell2mat(data_set10(:,14)) cell2mat(data_set15(:,14)) cell2mat(data_set20(:,14)) cell2mat(data_set25(:,14))];

inc_gain    = [cell2mat(data_set0(2,19)) cell2mat(data_set5(2,19)) cell2mat(data_set10(2,19)) cell2mat(data_set15(2,19)) cell2mat(data_set20(2,19)) cell2mat(data_set25(2,19))];
inc_control = [cell2mat(data_set_base(2,20)) ];
[inc_gain' inc_control'];
perc_gain = 100*(inc_gain'./inc_control'-1);


educ_gain    = [cell2mat(data_set0(2,34)) cell2mat(data_set5(2,34)) cell2mat(data_set10(2,34)) cell2mat(data_set15(2,34)) cell2mat(data_set20(2,34)) cell2mat(data_set25(2,34))]
educ_control = [cell2mat(data_set_base(2,35)) ]

eu_1p     = [cell2mat(data_set_base(2,16)), cell2mat(data_set0(2,16)), cell2mat(data_set5(2,16)),cell2mat(data_set10(2,16)),cell2mat(data_set15(2,16)), cell2mat(data_set20(2,16)), cell2mat(data_set25(2,16))]'
eu_change =-(eu_1p./eu_1p(1)-1)
eu_diff   = eu_1p-eu_1p(1)

spillover = [cell2mat(data_set_base(2,3)); cell2mat(data_set0(2,3)); cell2mat(data_set5(2,3)); cell2mat(data_set10(2,3)); cell2mat(data_set15(2,3)); cell2mat(data_set20(2,3)); cell2mat(data_set25(2,3))];
r_s = [cell2mat(data_set_base(2,5)); cell2mat(data_set0(2,5)); cell2mat(data_set5(2,5)); cell2mat(data_set10(2,5)); cell2mat(data_set15(2,5)); cell2mat(data_set20(2,5)); cell2mat(data_set25(2,5))];
n_s = [cell2mat(data_set_base(2,4)); cell2mat(data_set0(2,4)); cell2mat(data_set5(2,4)); cell2mat(data_set10(2,4)); cell2mat(data_set15(2,4)); cell2mat(data_set20(2,4)); cell2mat(data_set25(2,4))];

gini = [cell2mat(data_set_base(3,2)) cell2mat(data_set0(3,2)) cell2mat(data_set5(3,2)) cell2mat(data_set10(3,2)) cell2mat(data_set15(3,2)) cell2mat(data_set20(3,2)) cell2mat(data_set25(3,2))];
dissim = [cell2mat(data_set_base(2,1)) cell2mat(data_set0(2,1)) cell2mat(data_set5(2,1)) cell2mat(data_set10(2,1)) cell2mat(data_set15(2,1)) cell2mat(data_set20(2,1)) cell2mat(data_set25(2,1))];

% Partial Equilibrium
pe_data_set5  = counter_fact_mto(x_s, [], 0.05,1); % 5%
pe_data_set10 = counter_fact_mto(x_s, [], 0.10,1); % 10%
pe_data_set15 = counter_fact_mto(x_s, [], 0.15,1); % 15%
pe_data_set20 = counter_fact_mto(x_s, [], 0.20,1); % 20%
pe_data_set25 = counter_fact_mto(x_s, [], 0.25,1); % 25%
pe_data_set0  = counter_fact_mto(x_s, [], 0.0015,1); % MTO

inc_gain_pe = [cell2mat(pe_data_set0(2,19)) cell2mat(pe_data_set5(2,19)) cell2mat(pe_data_set10(2,19)) cell2mat(pe_data_set15(2,19)) cell2mat(pe_data_set20(2,19)) cell2mat(pe_data_set25(2,19))]
perc_gain_pe = 100*(inc_gain_pe'./inc_control'-1);

eu_1p_pe = [cell2mat(data_set_base(2,16)), cell2mat(pe_data_set0(2,16)), cell2mat(pe_data_set5(2,16)),cell2mat(pe_data_set10(2,16)),cell2mat(pe_data_set15(2,16)), cell2mat(pe_data_set20(2,16)), cell2mat(pe_data_set25(2,16))]'
eu_change_pe = -(eu_1p_pe./eu_1p_pe(1)-1)
eu_diff_pe = eu_1p_pe-eu_1p_pe(1)

%% Figures
scale = [0.15 5 10 15 20 25];

fig=figure;
hold on
plot(scale, perc_gain,'- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, perc_gain_pe,'-- .','LineWidth',1.5, 'MarkerSize',20)
lgd=legend({'General Equilibrium','Partial Equilibrium'}, 'Location','southoutside','FontSize',12,'Box','off');
lgd.NumColumns = 2;
set(gca,'XGrid','off','YGrid','on')
ax = gca;
ax.FontSize = 12; 
xtickformat('percentage')
ytickformat('percentage')
ylim([0 25])
yticks([0 5 10 15 20 25])
xlabel('Scale')
ylabel('Children % Income Gain')
hold off
set(fig,'Units','Inches');
pos = get(fig,'Position');
set(fig,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)])
print(fig,'figures/fig-incomegain.pdf','-dpdf','-r0')

fig=figure;
hold on
plot(scale, educ_gain,'- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, educ_control,'-- .','LineWidth',1.5, 'MarkerSize',20)
lgd=legend({'Treated','Control'}, 'Location','southoutside','FontSize',12,'Box','off');
lgd.NumColumns = 2;
set(gca,'XGrid','off','YGrid','on')
ax = gca;
ax.FontSize = 12; 
xtickformat('percentage')
yticks([0 0.05 0.10 0.15 0.20])
ylim([0 0.2])
xlabel('Scale')
ylabel('Education')
hold off
set(fig,'Units','Inches');
pos = get(fig,'Position');
set(fig,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)])
print(fig,'figures/fig-education.pdf','-dpdf','-r0')

scale = [0 0.15 5 10 15 20 25];

fig=figure;
hold on
plot(scale, spillover(:,1),'- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, spillover(:,2),'-- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, spillover(:,3),': .','LineWidth',1.5, 'MarkerSize',20,'color',[118,173,72]/255)
lgd=legend({'A','B','C'}, 'Location','southoutside','FontSize',12,'Box','off');
lgd.NumColumns = 3;
set(gca,'XGrid','off','YGrid','on')
ax = gca;
ax.FontSize = 12; 
xtickformat('percentage')
xlabel('Scale')
ylabel('Spillovers')
hold off
set(fig,'Units','Inches');
pos = get(fig,'Position');
set(fig,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)])
print(fig,'figures/fig-spillover-large.pdf','-dpdf','-r0')

fig=figure;
hold on
plot(scale, r_s(:,1),'- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, r_s(:,2),'-- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, r_s(:,3),': .','LineWidth',1.5, 'MarkerSize',20,'color',[118,173,72]/255)
lgd=legend({'A','B','C'}, 'Location','southoutside','FontSize',12,'Box','off');
lgd.NumColumns = 3;
set(gca,'XGrid','off','YGrid','on')
ax = gca;
ax.FontSize = 12; 
xtickformat('percentage')
xlabel('Scale')
ylabel('Rental Rates')
ylim([0.2 0.7])
hold off
set(fig,'Units','Inches');
pos = get(fig,'Position');
set(fig,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)])
print(fig,'figures/fig-rental-large.pdf','-dpdf','-r0')

fig=figure;
hold on
plot(scale, n_s(:,1),'- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, n_s(:,2),'-- .','LineWidth',1.5, 'MarkerSize',20)
plot(scale, n_s(:,3),': .','LineWidth',1.5, 'MarkerSize',20,'color',[118,173,72]/255)
lgd=legend({'A','B','C'}, 'Location','southoutside','FontSize',12,'Box','off');
lgd.NumColumns = 3;
set(gca,'XGrid','off','YGrid','on')
ax = gca;
ax.FontSize = 12; 
xtickformat('percentage')
xlabel('Scale')
ylabel('Neighborhood Size')
hold off
set(fig,'Units','Inches');
pos = get(fig,'Position');
set(fig,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)])
print(fig,'figures/fig-size-large.pdf','-dpdf','-r0')

